import matplotlib.pyplot as plt
import numpy as np

class IMCPID:
    def __init__(self, K=1.0, tau=50.0, theta=5.0, lam=40.0, dt=1.0):
        """
        IMC-PID控制器参数
        :param K: 系统增益
        :param tau: 系统时间常数
        :param theta: 系统时延
        :param lam: IMC滤波参数
        :param dt: 控制周期
        """
        # 计算 IMC-PID 参数
        self.Kp = (2 * tau + theta) / (K * (2 * lam + theta))
        self.Ti = tau + theta / 2
        self.Td = (tau * theta) / (2 * tau + theta)
        self.dt = dt

        # 控制器状态
        self.integral = 0.0
        self.prev_error = 0.0
        self.prev_measurement = 0.0

    def compute(self, setpoint, measurement):
        error = setpoint - measurement

        # 积分项 - 带抗饱和和积分限幅
        self.integral += error * self.dt
        self.integral = np.clip(self.integral, -5.0, 5.0)

        # 微分项 - 使用测量值微分减少设定值突变的影响
        derivative = (self.prev_measurement - measurement) / self.dt

        # IMC-PID 输出
        output = (self.Kp * error + 
                 (self.Kp / self.Ti) * self.integral - 
                 self.Kp * self.Td * derivative)
        
        self.prev_error = error
        self.prev_measurement = measurement
        return output


class HVACController:
    def __init__(self, season="winter"):
        self.season = season.lower()
        
        # 根据季节调整PID参数
        if self.season == "winter":
            self.pid = IMCPID(K=0.8, tau=60.0, theta=5.0, lam=40.0, dt=1.0)
        else:  # summer
            self.pid = IMCPID(K=1.2, tau=40.0, theta=3.0, lam=30.0, dt=1.0)
        
        # 控制参数
        self.min_on_time = 300  # 最小开启时间 (秒)
        self.min_off_time = 300  # 最小关闭时间 (秒)
        self.hysteresis = 0.3  # 迟滞带 (°C)
        self.emergency_shutdown_threshold = 2.0  # 紧急关闭阈值(°C)
        
        # 状态变量
        self.last_switch_time = 0
        self.current_state = "OFF"
        self.time_elapsed = 0

    def control(self, setTemp, roomTemp):
        pid_output = self.pid.compute(setTemp, roomTemp)
        self.time_elapsed += 1
        
        # 紧急关闭检查
        if (self.season == "winter" and roomTemp > setTemp + self.emergency_shutdown_threshold) or \
           (self.season == "summer" and roomTemp < setTemp - self.emergency_shutdown_threshold):
            self.current_state = "OFF"
            self.last_switch_time = self.time_elapsed
            return f"HVAC OFF (EMERGENCY)", pid_output
        
        # 根据季节确定控制方向
        temp_diff = setTemp - roomTemp
        if self.season == "winter":
            need_operation = temp_diff > self.hysteresis  # 需要加热
            stop_operation = temp_diff < -self.hysteresis
        else:  # summer
            need_operation = temp_diff < -self.hysteresis  # 需要制冷
            stop_operation = temp_diff > self.hysteresis
        
        # 检查最小运行/停机时间
        time_since_switch = self.time_elapsed - self.last_switch_time
        
        if self.current_state == "ON":
            if stop_operation and time_since_switch >= self.min_on_time:
                self.current_state = "OFF"
                self.last_switch_time = self.time_elapsed
        else:
            if need_operation and time_since_switch >= self.min_off_time:
                self.current_state = "ON"
                self.last_switch_time = self.time_elapsed
        
        return f"HVAC {self.current_state}", pid_output


class HouseModel:
    def __init__(self, initial_temp=18.0, season="winter"):
        self.temp = initial_temp
        self.season = season.lower()
        self.window_open = False  # 窗户状态
        self.window_open_start = 0
        self.window_open_duration = 0
        
        # 根据季节设置室外温度
        if self.season == "winter":
            self.outdoor_temp = 10.0  # 冬季室外温度
            self.heat_loss_coeff = 0.01  # 热损失系数
            self.heating_power = 0.5  # 加热功率
            self.cooling_power = 0.0  # 夏季制冷功率
        else:  # summer
            self.outdoor_temp = 32.0  # 夏季室外温度
            self.heat_loss_coeff = 0.015  # 夏季热损失更大(热空气进入)
            self.heating_power = 0.0  # 冬季加热功率
            self.cooling_power = 0.6  # 制冷功率
        
        self.thermal_mass = 100.0  # 热质量
        self.window_open_heat_loss = 0.1  # 开窗时的热损失系数
    
    def trigger_window_open(self, time, duration=300):
        """触发开窗事件"""
        self.window_open = True
        self.window_open_start = time
        self.window_open_duration = duration
    
    def update(self, hvac_on, time_step):
        # 检查开窗事件是否结束
        if self.window_open and (time_step - self.window_open_start >= self.window_open_duration):
            self.window_open = False
        
        # 热交换 - 与室内外温差成正比
        temp_diff = self.temp - self.outdoor_temp
        
        # 根据窗户状态调整热损失系数
        current_heat_loss_coeff = self.window_open_heat_loss if self.window_open else self.heat_loss_coeff
        heat_loss = current_heat_loss_coeff * temp_diff
        
        # HVAC效果
        if self.season == "winter":
            hvac_effect = self.heating_power if hvac_on else 0.0
        else:  # summer
            hvac_effect = -self.cooling_power if hvac_on else 0.0  # 制冷为负
        
        # 温度变化考虑热质量
        delta_temp = (hvac_effect - heat_loss) / self.thermal_mass
        self.temp += delta_temp
        
        # 添加微小随机扰动
        self.temp += np.random.normal(0, 0.005)
        
        # 确保温度不会超出物理极限
        self.temp = np.clip(self.temp, 0, 50)
        
        return self.temp, self.window_open


def simulate_hvac_with_dynamic_setpoints(season="winter", initial_temp=18.0, steps=10000):
    controller = HVACController(season=season)
    
    if season == "winter":
        house = HouseModel(initial_temp=initial_temp, season="winter")
        # 冬季设定温度变化序列: 22°C → 24°C → 28°C
        setpoints = np.ones(steps) * 22.0
        setpoints[steps//3:2*steps//3] = 24.0  # 中间1/3时间设为24°C
        setpoints[2*steps//3:] = 28.0  # 最后1/3时间设为28°C
    else:
        house = HouseModel(initial_temp=initial_temp+10, season="summer")
        # 夏季设定温度变化序列: 24°C → 22°C → 18°C
        setpoints = np.ones(steps) * 24.0
        setpoints[steps//3:2*steps//3] = 22.0  # 中间1/3时间设为22°C
        setpoints[2*steps//3:] = 18.0  # 最后1/3时间设为18°C
    
    temps = []
    actions = []
    outputs = []
    setpoint_history = []
    outdoor_temps = []
    window_events = []  # 记录窗户状态
    
    # 预定义开窗事件 (时间步长, 持续时间)
    window_events_def = [
        (steps//6, 300),     # 第一次开窗
        (steps//3 + 500, 600),  # 第二次开窗
        (2*steps//3 + 1000, 400) # 第三次开窗
    ]
    
    for t in range(steps):
        # 检查是否需要触发开窗事件
        for event_time, duration in window_events_def:
            if t == event_time:
                house.trigger_window_open(t, duration)
        
        # 模拟复杂的室外温度波动
        if season == "winter":
            # 冬季：基础温度 + 昼夜波动 + 随机天气变化
            base_temp = 10
            day_night = 5 * np.sin(t/1000)  # 昼夜波动
            weather_change = 3 * np.sin(t/3000)  # 缓慢的天气变化
            random_fluctuation = np.random.normal(0, 0.5)  # 随机波动
            house.outdoor_temp = base_temp + day_night + weather_change + random_fluctuation
        else:
            # 夏季：基础温度 + 昼夜波动 + 热浪 + 随机波动
            base_temp = 32
            day_night = 8 * np.sin(t/1000)  # 昼夜波动更大
            heat_wave = 5 * np.sin(t/5000)  # 热浪效应
            random_fluctuation = np.random.normal(0, 1.0)  # 更大的随机波动
            house.outdoor_temp = base_temp + day_night + heat_wave + random_fluctuation
        
        current_setpoint = setpoints[t]
        action, pid_output = controller.control(current_setpoint, house.temp)
        roomTemp, window_status = house.update("ON" in action, t)
        
        temps.append(roomTemp)
        actions.append(1 if "ON" in action else 0)
        outputs.append(pid_output)
        setpoint_history.append(current_setpoint)
        outdoor_temps.append(house.outdoor_temp)
        window_events.append(1 if window_status else 0)
    
    return temps, actions, outputs, setpoint_history, outdoor_temps, window_events


def plot_dynamic_setpoint_results(season, results):
    temps, actions, outputs, setpoints, outdoor_temps, window_events = results
    
    plt.figure(figsize=(14, 12))
    
    # 温度曲线
    plt.subplot(4, 1, 1)
    plt.plot(temps, label="Room Temperature", color="tab:blue", linewidth=1.5)
    plt.plot(setpoints, label="Setpoint", linestyle="--", color="tab:red", linewidth=1.5)
    plt.plot(outdoor_temps, label="Outdoor Temperature", color="tab:green", alpha=0.7)
    plt.title(f"{season.capitalize()} HVAC with Dynamic Setpoints and Disturbances")
    plt.ylabel("Temperature (°C)")
    plt.legend(loc="best")
    plt.grid(True)
    
    # PID输出
    plt.subplot(4, 1, 2)
    plt.plot(outputs, label="PID Output", color="tab:orange")
    plt.ylabel("PID Output")
    plt.legend(loc="best")
    plt.grid(True)
    
    # HVAC状态
    plt.subplot(4, 1, 3)
    plt.step(range(len(actions)), actions, label="HVAC State (1=ON,0=OFF)", 
             color="tab:green", where="post", linewidth=1.5)
    plt.ylabel("Control Signal")
    plt.legend(loc="best")
    plt.grid(True)
    
    # 窗户状态和扰动事件
    plt.subplot(4, 1, 4)
    plt.fill_between(range(len(window_events)), window_events, 
                     color="tab:purple", alpha=0.3, label="Window Open")
    plt.ylim(-0.1, 1.1)
    plt.xlabel("Time step")
    plt.ylabel("Disturbances")
    plt.legend(loc="best")
    plt.grid(True)
    
    plt.tight_layout()
    plt.show()


# 冬季动态设定温度仿真
print("Simulating winter scenario with dynamic setpoints (22°C → 24°C → 28°C)...")
winter_results = simulate_hvac_with_dynamic_setpoints(
    season="winter",
    initial_temp=18.0,
    steps=9000  # 3000步每个阶段
)

# 夏季动态设定温度仿真
print("Simulating summer scenario with dynamic setpoints (24°C → 22°C → 18°C)...")
summer_results = simulate_hvac_with_dynamic_setpoints(
    season="summer",
    initial_temp=28.0,
    steps=9000  # 3000步每个阶段
)

# 绘制结果
plot_dynamic_setpoint_results("winter", winter_results)
plot_dynamic_setpoint_results("summer", summer_results)